System and Method for Detecting Tagged Material Using Alpha Matting

ABSTRACT

A method for computer-aided object classification, soft segmentation and layer extraction in computed tomographic colonography includes providing a contrast enhanced computed tomography (CT) digital image of the colon, finding a foreground region of voxels with an intensity higher than a pre-defined foreground threshold, creating a 3D trimap of the colon where the image is segmented into the foreground region, a background region, and an unknown region between the foreground and background, starting from the background, extracting successive layers of the unknown region until the foreground region is reached, and classifying each extracted layer as background or foreground, and generating a foreground matte, a background matte, and an alpha matte, where alpha indicates a mixing ration of foreground and background voxels.

CROSS REFERENCE TO RELATED UNITED STATES APPLICATIONS

This application claims priority from “Alpha Matting Technique for theDetection of Tagged Material”, U.S. Provisional Application No.60/978,896 of Lee, et al., filed Oct. 10, 2007, the contents of whichare herein incorporated by reference in their entirety.

TECHNICAL FIELD

This disclosure is directed to object classification in computedtomographic colonography.

DISCUSSION OF THE RELATED ART

Colorectal cancer is the second leading cause of cancer-related death inthe western world. Early detection of precancerous polyps throughcolorectal screening can largely prevent colon cancer by enablingremoval of the polyps before they turn malignant. However, traditionaloptical colonoscopy is highly uncomfortable which hinders the popularityof routine colorectal screening. Virtual colonoscopy, on the other hand,uses minimally invasive computer tomography (CT) scanning which reducesthe cost and risk and still yields a high diagnostic performance todetect clinically significant polyps as compared to traditionalcolonoscopy. However, since residual materials inside the colon have thesame intensity as the colonic wall in a CT scan, they could be falselyinterpreted as a part of the colon which can hinder the detection ofpolyps. Some polyps may remain undetected and some residual material maybe misinterpreted as polyps. Therefore, proper bowel preparation knownas colon cleansing is carried out. This usually includes cleansing thecolon with large amount of fluid or administering laxatives to inducebowel movements. The process is uncomfortable and is a major bottleneckto the adoption of CT colon screening. Several groups of researchershave developed different fluid-based diets for bowel preparation. Afterthe diet, the colon of the patient contains mostly residual fluid. Withthe ingestion of contrast agents, the intensity of the residual fluidwill be enhanced. By detecting and removing these high intensityresiduals from the image, the colon image is cleared from residuals andpolyp detection can then be carried out. However, patients still need toundergo several days of a fluidic diet which is still uncomfortable.Because of this, the colon contains not only residual fluids but alsosolid residuals attached along the colonic wall and folds. These solidresiduals hinder the detection of polyps and can result in falsepositives. However, due to the limited resolution and contrast of CT andthe varying absorption rate of the contrast agents, the enhancedresidual fluids do not have clear boundaries due to partial voluming,which hinders the segmentation and removal of these residuals. Theboundary between air and residuals has a varying thickness and theintensity of the boundary is the same as soft tissue. In addition, thesoft tissue near the boundary between residuals and tissue is enhancedand is much higher than the typical intensity of soft tissue.

To remove these solid tagged materials, there are several challenges:

1. Partial voluming and pseudo-enhancement effects due to the highintensity of tagged materials.

2. The intensity of tagged materials can range from 1300HU up to 4000HU,with different effects on partial volume and pseudo-enhancement due todifferent intensities.

3. Due to different absorption rates of the contrast agents, theintensity of tagged material drops along the colon, with differentintensities in different locations of the colon.

4. Solid and liquid tagged materials can form a thin layer along thecolonic wall. The severe partial voluming can hinder their removal.

Because of the high intensity of tagged materials, they can create anaura around themselves that causes an intensity glow on the neighborhoodelements, which is called pseudo-enhancement. Traditional methods ofusing thresholding to remove tags from air and tissue will fail in thiscase. The thickness and intensity of the aura varies greatly, from 0 to5 voxels, due to inhomogeneous absorption of the contrast agent by thetagged material. The intensity of the tagged material can range from1100HU to 5000HU. Extracting and removing unwanted objects or occlusionfrom 2D image and 3D scenes has been an important subject in computervision, and is especially important in medical imaging, where occlusionsmay hinder the detection of important objects valuable for diseasedetection in 3D tomography data. To accurately extract an object from a3D dataset, one should to take into account the “intensity blending”across object boundaries to recover the pure intensity of thecorresponding object and background in the ambiguous blended regions.

To understand the challenges involved in the removal of colon solidtagged materials, referred to herein below as tags, consider thesimplest and traditional thresholding method and see how it fails.

FIG. 1 illustrates the effects of performing thresholding on an image.If one sets a high threshold, a hole 11 will be formed in the tag, witha thin line remaining. If one sets a low threshold, some tissue voxelswill be removed.

Current research on image segmentation in the presence of occlusions hasnot addressed the issue of opacity in the presence of partial volumes.However, it is the blending ratio of the partial volume that is ofinterest, not the recovery of the occluded part, one doe not need toinpaint the lost region using level lines or intensity information.

By finding the opacity in the presence of pixel color blending, mattingis a powerful tool for natural foreground extraction. Matting typicallyrequires a user to specify a trimap which indicates a definiteforeground where α equals 1, a definite background where α equals 0, andan unknown region where α is to be determined. Information from knownbackground and foreground is utilized to estimate the α in the unknownregion. Matting can obtain a soft segmentation based on the α values. Itis particularly useful when images have severe partial volumes andlimited resolution. However, currently there are no matting algorithmsbeing applied to medical tomographic images, partly due to the fact thatmatting was originally designed for extracting a foreground such as ahuman face to “matte” it to another background. But in fact, matting andtrimap creation provide powerful insights for handling partial volumingand pseudo-enhancement tasks in 3D colon cleansing.

Tissue classification is also an active area of research in MRI and CT.Current methods focus on multiple tissue types which do not have cleargeometric relationship, such as in the colon where soft tissue, air andtag locations can be easily modeled.

Most of the previous works on colon cleansing has focused on taggedliquid materials after the liquid diet. Because the tagged materials areliquid, they show a horizontal line parallel to the bed where thepatient is lying on. Above the liquid level the voxels are definitelybelonged to air. The most challenging task is to locate that boundary ofliquid and air to remove the liquid.

SUMMARY OF THE INVENTION

Exemplary embodiments of the invention as described herein generallyinclude methods and systems for fully-automatic 3D computer-aided objectclassification, soft segmentation and layer extraction in the newestcomputed tomographic digital colon fecal-tagging protocol for earlydetection of colorectal cancer, by unifying and generalizing thetechniques in alpha matting, morphological operations and featureaggregation. A method according to an embodiment of the inventionobtains an alpha value which indicates the blending of foreground andbackground, and accurately classifies the background to resolve partialvoluming and pseudo-enhancement effects, and extracts and removesresidual solid materials for patients who do not undergo any bowelpreparation diet. A method according to an embodiment of the inventionincludes a new method called alphamat detagging, which extends mattingto 3D tomographic images, fully automatic trimap generation in medicalimaging through anatomical and geometric priors without user input,automatic classification of tissue types in the presence of partialvolume and pseudo-enhancement which hinders the automated detection ofpolyps in fecal-tagging colon computed tomography. A confidence measureis introduced which assigns a probability to each voxel of the datasetincluding the regions where foreground objects are fractionally blendedwith the background objects. Experimental results are presented whichshow that a method according to an embodiment of the invention canextract ambiguous boundaries such as foreground objects with shapes of“thin cloth” or “small balls” where classical methods are notapplicable. The results indicate that laxatives and a prescreening fluiddiet may no longer be needed for colon cancer CT screening. An algorithmaccording to an embodiment of the invention can also automaticallydetect T-junctions on-the-fly.

According to an aspect of the invention, there is provided a method forcomputer-aided object classification, soft segmentation and layerextraction in computed tomographic colonography, including providing acontrast enhanced computed tomography (CT) digital image of the colon,the image comprising a plurality of intensities associated with a3-dimensional (3D) grid of voxels, finding a foreground region of voxelswith an intensity higher than a pre-defined foreground threshold,creating a 3D trimap of the colon where the image is segmented into theforeground region, a background region, and an unknown region betweenthe foreground and background, starting from the background, extractingsuccessive layers of the unknown region until the foreground region isreached, and classifying each extracted layer as background orforeground, and generating a foreground matte, a background matte, andan alpha matte, where alpha indicates a mixing ration of foreground andbackground voxels.

According to a further aspect of the invention, the background regionvoxels represent either regions of air in the colon, or regions of knownsoft tissue in the colon.

According to a further aspect of the invention, creating a 3D trimap ofthe colon comprises morphologically dilating the 3D foreground regionone voxel at a time until background region voxels are reached, where aregion through which the foreground was diluted is a region of unknowntype, and categorizing the background region voxels as either air orsoft tissue.

According to a further aspect of the invention, classifying eachextracted layer as background or foreground comprises forming a 3Dsphere around a candidate layer voxel and neighboring known backgroundvoxels, and analyzing closest known background voxels to determinewhether the candidate layer voxel is a background air voxel or abackground soft tissue voxel.

According to a further aspect of the invention, the method includescalculating the mixing ratio alpha for each voxel in the unknown regionfrom alpha=(C−B)/(F−B), where C is the intensity of the voxel, B is abackground intensity, and F is a foreground intensity.

According to a further aspect of the invention, analyzing closest knownbackground voxels comprises finding values of f and B that maximize asum of log-likelihoods L(C|B, F, α)+L(F)+L(B)+L(α), where

${{L\left( {{CF},B,\alpha} \right)} = {- \frac{{{C - {\alpha \; F} - {\left( {1 - \alpha} \right)B}}}^{2}}{2\sigma_{C}^{2}}}},{{L(F)} = {- \frac{{{F - \overset{\_}{F}}}^{2}}{2\sigma_{F}^{2}}}},{{{and}\mspace{14mu} {L(B)}} = {- \frac{{{B - \overset{\_}{B}}}^{2}}{2\sigma_{B}^{2}}}},$

where σ_(C) ^(r), σ_(F) ², σ_(B) ² are respective standard deviationsfor C, F, and B.

According to a further aspect of the invention, generating a foregroundmatte, a background matte, and an alpha matte, comprises forming a 3Dsphere about each candidate voxel in the region of unknown type,obtaining from the background region within the sphere a number of airvoxels, a number of tissue voxels, their respective distance from thecandidate voxel, and their variance, obtaining from the foregroundregion within the sphere a number of tagged voxels having an intensityabove a threshold for contrast enhancement, their respective distancefrom the candidate, and their variance, weighting the intensity of eachtissue voxel by its distance from the candidate voxel, weighting theintensity of each air voxel by its distance from the candidate voxel,weighting the intensity of each foreground voxel by its distance fromthe candidate voxel, and calculating alpha fromC=αF+∫_(Ω(air))H_(air)(1−α(i)B_(air)(i)+∫_(Ω(tissue))H_(tissue)(1−α(i))B_(tissue)(i),where C is the intensity of the candidate voxel, B_(air) is a backgroundair intensity, B_(tissue) is a background tissue intensity, F is aforeground intensity, Ω is a region of integration in either air voxelsor tissue voxels, i represents a voxel in the region of integration, andH is the Heaviside function having a value of 1 inside the region ofintegration.

According to a further aspect of the invention, a foreground weightingfunction is

${F_{candidate} = {{\exp \left( {{- \beta_{1}}{\sum d_{i}}} \right)}{\exp \left( {{- \beta_{2}}\text{/}N} \right)}\frac{\sum{F_{i}\text{/}d_{i}}}{\sum{1\text{/}d_{i}}}}},$

where F_(candidate) is the foreground weighted intensity at thecandidate voxel, F_(i) is the intensity value of the image at voxel iwithin the 3D sphere, d_(i) is the distance to the candidate voxel, β₁and β₂ are pre-defined parameters, a background weighting function is

${B_{candidate} = {{\exp \left( {{- \beta_{1}}{\sum d_{i}}} \right)}{\exp \left( {{- \beta_{2}}\text{/}N} \right)}\frac{\sum{B_{i}\text{/}d_{i}}}{\sum{1\text{/}d_{i}}}}},$

where B_(candidate) is the background weighted intensity at the voxellocation and B_(i) indicates the intensity value of the known backgroundwithin the 3D sphere.

According to a further aspect of the invention, a candidate voxel isidentified as a T-junction where the number of tissue voxels, airvoxels, and a number of tagged voxels are greater than zero inside asphere of a predetermined small radius.

According to a further aspect of the invention, for tagged voxels layingon a haustral fold, the method includes dilating the tagged voxels by anumber of voxels less than an original size of the morphologicaldilation, removing the tagged voxels and dilated voxels from a voxelwindow whose side is greater than the original size of the morphologicaldilation, performing region growing to find connected voxels,identifying the connected voxels, and performing matting, where there isno type change of any voxel that has been identified as a connectedvoxel.

According to another aspect of the invention, there is provided aprogram storage device readable by a computer, tangibly embodying aprogram of instructions executable by the computer to perform the methodsteps for computer-aided object classification, soft segmentation andlayer extraction in computed tomographic colonography.

BRIEF DESCRIPTION OF THE DRAWINGS

FIG. 1 illustrates how thresholding can create holes in an object,according to an embodiment of the invention.

FIG. 2 illustrates partial volume and pseudo enhancement effects,according to an embodiment of the invention.

FIG. 3 illustrates tomographic matting with a trimap, according to anembodiment of the invention.

FIG. 4. is a flowchart of a tomographic matting method for colondetagging, according to an embodiment of the invention.

FIG. 5 illustrates an ‘onion peeling’ that marches inwards from a knownbackground, according to an embodiment of the invention.

FIG. 6 depicts a 2D example of onion peeling, according to an embodimentof the invention.

FIG. 7 illustrates how a 3D sphere captures information from backgroundand foreground relative to a candidate location, according to anembodiment of the invention.

FIGS. 8( a)-(d) illustrate a tag on a fold, according to an embodimentof the invention.

FIGS. 9( a)-(d) show the results of classification and tag removal froma thin cloth shape, according to an embodiment of the invention.

FIGS. 10( a)-(d) illustrate a real tissue voxel appearing to be partialvolume between air and tag, according to an embodiment of the invention.

FIGS. 11( a)-(d) show the effect of a visible tag below the thresholdintensity, according to an embodiment of the invention.

FIGS. 12( a)-(c) show a case where there is a hole between the tag andthe tissue, according to an embodiment of the invention.

FIG. 13 is a block diagram of an exemplary computer system forimplementing a method for detecting tagged material using alpha mattingin computed tomographic colonography, according to an embodiment of theinvention.

DETAILED DESCRIPTION OF EXEMPLARY EMBODIMENTS

Exemplary embodiments of the invention as described herein generallyinclude systems and methods for detecting tagged material using alphamatting in computed tomographic colonography. Accordingly, while theinvention is susceptible to various modifications and alternative forms,specific embodiments thereof are shown by way of example in the drawingsand will herein be described in detail. It should be understood,however, that there is no intent to limit the invention to theparticular forms disclosed, but on the contrary, the invention is tocover all modifications, equivalents, and alternatives falling withinthe spirit and scope of the invention.

As used herein, the term “image” refers to multi-dimensional datacomposed of discrete image elements (e.g., pixels for 2-D images andvoxels for 3-D images). The image may be, for example, a medical imageof a subject collected by computer tomography magnetic resonanceimaging, ultrasound, or any other medical imaging system known to one ofskill in the art. The image may also be provided from non-medicalcontexts, such as, for example remote sensing systems, electronmicroscopy, etc. Although an image can be thought of as a function fromR³ to R, the methods of the inventions are not limited to such images,and can be applied to images of any dimension, e.g., a 2-D picture or a3-D volume. For a 2- or 3-dimensional image, the domain of the image istypically a 2- or 3-dimensional rectangular array, wherein each pixel orvoxel can be addressed with reference to a set of 2 or 3 mutuallyorthogonal axes. The terms “digital” and “digitized” as used herein willrefer to images or volumes, as appropriate, in a digital or digitizedformat acquired via a digital acquisition system or via conversion froman analog image.

Overview

Matting has successfully evolved in computer graphics with results thatare well-validated by natural images, high dynamic image ground truth,and synthetic images. Alpha matting is used to estimate the ambiguousblended region at the boundary of object and background. Since colonsolid tag can have any size and shape, and be located anywhere along thecolonic wall and folds, machine learning, hierarchical or multiscalemethods are not suitable here. Also, since the tag has varyingintensity, size and shape, and can light up nearby voxels, hardsegmentation with a definite boundary is not a good approach.

Instead, an alpha channel is used which captures the idea that theintensity value in an image voxel can be a result of mixing ofintensities from multiple objects, resulting in a need for softsegmentation. The color blending is due to limited resolution of CT. Inthis framework, one looks not only at the neighborhood intensity, but ineach cluster of neighborhoods, one looks at a rich set of featuresincluding distance, number of neighborhood voxels, their gradients,variance, etc. A 3D sphere is created to handle partial voluming basedon the rich set of features of the neighborhood voxels. Since thebackground information is propagated layer by layer from the regionwhich is definitely not contaminated with pseudo-enhancement, the taskof pseudo-enhancement is addressed at the same time.

No seed point is needed to start the segmentation because region growingis not needed as initial information in colon. No learning is needed toanalyze the interface between tagged material/air, air/tissue,tissue/tagged material as everything can be determined by knownbackground propagation. Unknown regions are defined by two criteria.First, soft tissue will be mistakenly detected at the mixed area oftagged material and air. Secondly, intensity of real soft tissue will beincreased at the interface of tissue and tagged material. A confidencelevel is assigned to each voxel at the mixed boundary to indicate thecontribution of foreground intensity and background intensity to giverise to the intensity of the mixed voxel. The task is inherentlyunderconstrained because for each observation C, one needs to find itscorresponding foreground F, background B, and ratio. In computergraphics with color images having complicated backgrounds andforegrounds, this task is ill-posed. However, a set of prior informationcan be used. For example, one knows that tagged material is alwaysinside air, on the colonic wall or folds. FIG. 2 illustrates partialvolume (PV) and pseudo-enhancement (PE) effects in an image. On the leftis a simulated model showing a colon wall 21, tag 22, and PV and PE 23.On the right is real data with colon wall 21, tag 22, and PV and PE 23.The grey area 23 around the white tag is the aura, which creates thepartial voluming and pseudo-enhancing.

A method according to an embodiment of the invention can assignprobabilities to indicate the blending ratio of solid residual materialsand colonic wall, where partial voluming and pseudo enhancement aresevere, and provide a confidence map of residual materials to theexisting colon CT software to further improve the accuracy of polypdetection. A method according to an embodiment of the invention is afirst application of alpha matting to colon cancer detection, and firstapplication of assigning a confidence probability to a foreground toaddress the solid tag issue for a protocol where patients do not need toundergo fluid diet prior to a scan.

Tomographic Matting Model

A composition equation according to an embodiment of the invention is:

C=αF+(1−α)B,  (1)

where C, F, and B are a voxel's composite, foreground and backgroundintensity, respectively, and α is a fractional opacity with a valuebetween 0 and 1 which is the mixing ratio of the foreground intensityand background intensity. The basic idea is to estimate intensity F, andB, and the ratio α from the observed intensity C.

There is one grayscale image matting equation (EQ. (1)) and threeunknowns. Given the intensity of C, a Bayesian matting method finds themost likely values of F, B, and α by maximizing over a sum oflog-likelihoods:

$\begin{matrix}{{\underset{{F.B},\alpha}{argmax}{P\left( {F,B,{\alpha C}} \right)}} = {\underset{F,B,\alpha}{argmax}\frac{P\left( {{CF},B,\alpha} \right){P(F)}{P(B)}{P(\alpha)}}{P(C)}}} \\{= {\underset{F,B,\alpha}{argmax}\begin{pmatrix}{{L\left( {{CB},F,\alpha} \right)} +} \\{{L(F)} + {L(B)} + {L(\alpha)}}\end{pmatrix}}}\end{matrix}$

where L( ) is the log-likelihood of P( ). The first term introduces astandard deviation σ_(C) and is used to model the Gaussian noise in theimage:

$\begin{matrix}{{L\left( {{CF},B,\alpha} \right)} = {- {\frac{{{C - {\alpha \; F} - {\left( {1 - \alpha} \right)B}}}^{2}}{2\sigma_{C}^{2}}.}}} & (2)\end{matrix}$

By capturing the foreground voxels within a circular radius, the voxelsare partitioned into several clusters, that is, into connected singleregions. For each cluster, which is a combination of colors in colormatting, a weighted mean color F and weighted covariance matrix Σ_(F)are calculated. The weight of a voxel is related to α_(i) ², which givesthe color of a more opaque voxel higher confidence. In a colon CT, amore opaque voxel has a higher chance of being tagged material. Aspatial Gaussian g_(i) can be used to stress the contribution of nearbypixels over more distant voxels. An exemplary, non-limiting value of σis 8 voxels. The combined weight is w_(i)=α_(i) ²g_(i):

$\begin{matrix}{{\overset{\_}{F} = {\frac{1}{W}{\sum\limits_{i \in N}{w_{i}F_{i}}}}},} & (3) \\{{\Sigma_{F} = {\frac{1}{W}{\sum\limits_{i \in N}{\left( {F_{i} - \overset{\_}{F}} \right)\left( {F_{i} - \overset{\_}{F}} \right)^{T}}}}},} & (4)\end{matrix}$

There is a detailed discussion of the weighted distance function below.The log likelihoods for the foreground L(F) can be modeled as:

$\begin{matrix}{{{L(F)} = {{- \frac{1}{2}}\left( {F - \overset{\_}{F}} \right)^{T}{\Sigma_{F}^{- 1}\left( {F - \overset{\_}{F}} \right)}}},} & (5)\end{matrix}$

where Σ_(F) ⁻¹ is the inverse covariance matrix. For grayscale image, itbecomes:

$\begin{matrix}{{L(F)} = {\frac{{{F - \overset{\_}{F}}}^{2}}{2\sigma_{F}^{2}}.}} & (6)\end{matrix}$

The same applies to estimating background:

$\begin{matrix}{{L(B)} = {- {\frac{{{B - \overset{\_}{B}}}^{2}}{2\sigma_{B}^{2}}.}}} & (7)\end{matrix}$

FIG. 3 illustrates tomographic matting with a trimap, illustrating thegeometric properties of foreground 31, background soft tissue 33 andbackground air 34, and unknown tissue 32 in a colon CT. This allowsbreaking the different intensity distributions of the clusters into 3types, thus allowing fast computations. By exploiting these geometricproperties, there is no need to compare and choose the closest clustersas candidates for foreground and background because one geometricallyforces a single cluster for background and foreground respectively.Without the computationally intensive comparison of differentcandidates, and without the need to project colors because of thegrayscale image, one can directly obtain α based on weighted mean offoreground and weighted mean of background:

$\begin{matrix}{{C = {{\alpha \; F} + {\left( {1 - \alpha} \right)B}}},} \\{= {{{\alpha (i)}{F(i)}} + {\int_{\Omega {({air})}}^{\;}{{H_{{air}/{tag}}\left( {1 - {\alpha (i)}} \right)}{B_{air}(i)}}} +}} \\{{\int_{\Omega {({tissue})}}^{\;}{{H_{{tissue}/{tag}}\left( {1 - {\alpha (i)}} \right)}{{B_{tissue}(i)}.}}}}\end{matrix}$

where Ω indicates the region with the respective blended background, andH is the Heaviside function, which is 1 in the region indicated and 0elsewhere.

Two challenging aspects in previous matting methods are creating trimapsthrough user interaction, and clustering the foreground and backgroundinto oriented Gaussians with weighted mean and covariance matrices andselecting the contributions that maximize the equation.

FIG. 4 shows a flowchart of an algorithm according to an embodiment ofthe invention. An algorithm begins at step 41, when a 3D foreground iscreated by finding voxels with an intensity higher than a predefinedthreshold for tagged material. An exemplary, non-limiting thresholdvalue is 1300HU. Then, at step 42, a 3D trimap is created. First, amorphological dilation of a the 3D foreground region from step 41 isperformed to identify the unknown region, based on the geometricproperties of colon that the unknown region is always the “aura”/partialvolume surrounding the foreground. An exemplary, non-limiting dilationwindow is 7×7×7, with support in the middle, and a thickness of theunknown material is from 0 to 4 voxels. The foreground region is dilatedone voxel at a time to a width larger than the definite partial volumeregion. Prior information of the intensity values of tissue and/or airis used along with the gradient values to stop the dilation. So, if airor tissue voxels are encountered, then the dilation is stopped in thatdirection. The region outside the dilated region is definite background.After the foreground and unknown matte are formed, one can furtheridentify the known, background region into definite-tissue ordefinite-air by using prior knowledge of the air threshold and the softtissue threshold. After these steps, an accurate 3D trimap is formedthat segments the colon into a foreground, a background, and an unknownregion. At step 43, onion-peeling is performed by eroding the image onelayer at a time, propagating from the background side of the unknownregion to the foreground side. At each layer, the real tissue type isdetermined based on the previous layer using a Bayesian matting modelaccording to an embodiment of the invention described above, whicheliminates partial volume and pseudo-enhancement effects. The alpha andbackground are recovered from the foreground at step 44. At each layer,a 3D sphere of a predefined radius is used to calculate foreground andbackground statistics, and the weighted means and variances of clustersin the sphere are calculated. Then, alpha can be calculated and thebackground, foreground, and alpha matte can be generated.

By extending the Bayesian matting to 3D, a method according to anembodiment of the invention can determine foreground color F, backgroundcolor B, and opacity given the observed intensity C for each voxelwithin the unknown region of the 3D dataset. A method according to anembodiment of the invention can automatically create trimaps withoutuser-input, and uses an “onion peel” continuously sliding window tode-fine neighborhoods, propagating inward from the background regionsbut not from the known foreground, due to the ambiguity of tissue typenear the foreground, which can only be solved from known informationpropagated from background. An algorithm according to an embodiment ofthe invention utilizes nearby computed B, F and α through their Gaussiandistributions.

Automatic 3D Trimap Segmentation

Trimap determination in videos requires manual construction per frame,which is time consuming, and estimation of trimaps is sometimesinaccurate. But in medical tomographic images, there is a rich knowledgeof anatomy priors which can be applied to obtain an accurate trimap.Most layer extraction algorithms generate binary masks for classifyingpixels as foreground and background with a hard segmentation, instead ofusing opacity/alpha values to indicate the blending of colors andpartial volumes.

Based on the prior information of CT, where the fecal-tagged materialshave intensity higher than 1300HU, one can in an embodiment of theinvention set the threshold for foreground as 1300HU. Any voxel above1300HU is classified as definite Foreground:

$\begin{matrix}{{F\left( {i,j,k} \right)} = \left\{ \begin{matrix}1 & {{{{if}\mspace{14mu} {I\left( {i,j,k} \right)}} > {threshold}},} \\0 & {{otherwise}.}\end{matrix} \right.} & (8)\end{matrix}$

where I is the 3D CT image dataset, and F is the 3D definite foregroundmask. Note that this value is exemplary and non-limiting, and othervalues are within the scope of other embodiments of the invention.

Since there is prior information that F has higher intensity than B andF also has higher intensity than C, one can write

${\alpha = \frac{C - B}{F - B}},$

where B is either B_(air) or B_(tissue) depending on the backgroundpropagation. Because of the incorporation of medical priors, the taskreduces to a case where there is no need to search the best cluster forbackground. An onion peel method according to an embodiment of theinvention ensures that the background tissue type is always chosen basedon the nearest neighborhood of known background. The foregroundthresholding can be extended to more general case where other priorinformation about the objective foreground is known.

Adaptive Morphological Dilation

Due to the fact that partial volume and pseudo-enhancement around thetagged materials do not have constant thickness, one may need to computeredundant voxels that are unaffected by the tagged materials if onealways dilates by a constant number of voxels. However, when performingmorphological dilation, one can incorporate CT priors and geometricpriors to adaptively dilate to the unknown region which does not possessinformation regarding its background type. When dilating to a voxel thathas intensity lower than the air-threshold, the dilation is stopped andthat voxel is labeled as definite background and classified as air. Whendilating to a voxel that has near constant gradients in its x, y and zdirections and its intensity is in the soft tissue range, its truetissue is probability not affected by tagged materials.

Marching Inwards from Known Background

A method according to an embodiment of the invention marches from theknown background. The unknown region closest to the foreground in colonCT is ambiguous and it is challenging to determine a candidate type fromthat region. Instead of a sliding 2D window, a 3D sphere with apredefined radius is used as a window to gather neighborhoodinformation.

FIG. 5 illustrates an onion peeling mechanism according to an embodimentof the invention, where layers are extracted by marching inward from aknown background 51 to unknown tissue 52. A new candidate layer isinvestigated by a 3 by 3 by 3 kernel image erosion. At each candidatevoxel, the closest known background voxel is analyzed to determine thetype of candidate voxel, either tissue or air. FIG. 6 illustrates a 2Dexample of marching inwards from a known background, proceeding from theleft image to the right image.

3D Feature Aggregation and Scale Selection by Radius of Sphere

A method according to an embodiment of the invention uses a 3D sphere tocapture the neighborhoods of a candidate unknown voxel. 3D informationis more important than a slice by slice treatment because the underlyingstructures may be misleading in 2D but may be obvious in 3D. Anexemplary, non-limiting radius for colon cleansing is 4 voxels. The sizeof the radius indicates how much neighborhood information is desired.For a noisier dataset, a larger radius is useful to cancel noise. In acomplicated structure where tissue, air and tagged materials arepresent, a larger radius can analyze the best tissue classification,obtaining a coarser scale analysis. What can be obtained from the 3Dsphere neighborhood is a set of features, flexible for differentpurposes. In the colon cleansing case, one can obtain from the definitebackground region within the sphere radius the number of air voxels, thenumber of tissue voxels, their respective distance from the candidatevoxel, and their variance, which can be used for noise reduction. Fromthe definite foreground region within the sphere radius, one can obtainthe number of tagged voxels above the tag threshold, their respectivedistance from the candidate, and the variance of the tagged voxels.Other statistical and geometric measures can be obtained, based ondifferent applications.

FIG. 7 illustrates how a 3D sphere (drawn as 2D in the figure) capturesinformation from background and foreground relative to the candidatelocation. This figure shows how the sphere neighborhood can joinbackground 71, foreground 72, and unknown regions 73.

Classification and T-Junction

A 3D sphere approach according to an embodiment of the invention canalso obtain the location of a Tjunctions on the fly, while counting thenumber of background and foreground type voxels. When the number oftissue voxels, air voxels and tagged material voxels are greater thanzero but have a small radius, that candidate point is a T-junctionlocation.

Matting and Alpha Value

After the candidate type of each layer of an “onion peel” is determined,the 3D sphere neighborhood of each candidate is analyzed. If thecandidate type is tissue, the weighted intensity of tissue from the 3Dsphere neighborhood of definite background is calculated from theintensity of each tissue voxel in the background neighborhood weightedby its distance from the candidate. The weight doubles as the distancedoubles. The same applies to the case when candidate type is air. Thealpha value α determines how likely that a voxel belongs to taggedmaterial and how likely it is background. By tuning the segmentationbased on thresholding a particular value of alpha, there is a potentialthat the true background boundary can be recovered.

Weighted Distance Function

According to an embodiment of the invention, there is a weighteddistance measure which is proportional to the number of same-typebackground voxels in the 3D sphere, and their relative distance to thecandidate. A user can fine tune the weight. A linear blendingrelationship of foreground and background is assumed. The intensitydrop-off of foreground elements with respect to distance from thecandidate voxel should be related to the number of foreground voxels inthe 3D sphere, their relative distances to the candidate, and theirdistance-weighted mean intensity. A foreground weighting functionaccording to an embodiment of the invention is:

$\begin{matrix}{{F_{candidate} = {{\exp \left( {{- \beta_{1}}{\sum d_{i}}} \right)}{\exp \left( {{- \beta_{2}}/N} \right)}\frac{\sum{F_{i}/d_{i}}}{\sum{1/d_{i}}}}},} & (9)\end{matrix}$

where F_(candidate) indicates the foreground weighted intensity at thecandidate location, and F_(i) indicates the intensity value of the imageat location i within the 3D sphere. β₁ and β₂ are parameters that can beset by a user. Similarly, the background weighting function is:

$\begin{matrix}{{B_{candidate} = {{\exp \left( {{- \beta_{1}}{\sum d_{i}}} \right)}{\exp \left( {{- \beta_{2}}/N} \right)}\frac{\sum{B_{i}/d_{i}}}{\sum{1/d_{i}}}}},} & (10)\end{matrix}$

where B_(candidate) indicates the background weighted intensity at thecandidate location and B_(i) indicates the intensity value of the knownbackground within the 3D sphere which is propagated by onion peeling.

In one exemplary, non-limiting embodiment of the invention, β₁ and β₂were set to 0 for simplicity. The results will be more accurate if thesetwo parameters are fine-tuned.

More sophisticated information can be added to this relationship ifthere is more information regarding the foreground and background suchas a true profile between tagged material and soft tissue.

Soft Tissue, Air, and Partial Volumes

Assume a colon without tagged material. Inside the colon lumen, which iscomposed of soft tissue, there should only be air. Since it is desiredto remove the foreground and classify the restored background intoeither air or tissue, the resulting background will include only air andtissue and the partial volumes at the boundaries originally containingtagged materials will be gone. This may look unnatural, but for theidentification of true tissue from fake tissue to aid the polypdetection, visualization and cosmetic effects can be ignored. Theseboundaries can be smoothed or the partial volume effect can berecalculated based on the feature information aggregated by the 3Dsphere in those boundary candidate voxels.

Thresholded Foreground

Since the foreground region, which has an intensity higher than thethreshold, always has α equal to 1, this region will be considered to bepure air after the foreground is removed. However, in 3D tomography, thethresholded region sometimes covers tissue. If these regions areremoved, they could be replaced by tissue. However, it is not necessaryto do so, as one is more interested in the probability of tags in theunknown re-ion to avoid false positives in polyp detection.

Multiple Nearby Foreground Regions

If several foreground regions are nearby, the dilation process maycombine and close several unknown regions. The onion peeling may notreach some regions that are closed. For example, if two tags are 10voxels apart, then a 5 voxel dilation will completely close the gaps. Inthis situation, an onion peeling according to an embodiment of theinvention is still performed. If the tissue type is determined but theknown backgrounds are out of reach of the predefined radius of the 3Dsphere, the radius can be increased until it obtains some neighborhoodinformation.

Tagged Materials on Haustral Folds

One challenge in tag removal concerns tagged material laying on ahaustral fold. Since folds can appear as narrow as 4 mm in a CT slice, afold can have the same thickness and appearance as a partial volumebetween tagged materials and air. Since a matting algorithm according toan embodiment of the invention uses known neighborhood information, theknown neighborhood information of a haustral fold, which is mostly air,is not helpful in distinguishing partial volume from true tissue. Amatting algorithm will break the fold while removing the tag. However,according to an embodiment of the invention, instead of dilating thethresholded tag by the current number voxels, dilate it by a reducednumber of voxels, then subtract the thresholded tag and the dilated partfrom a voxel window whose side is greater than the original dilationsize. For example, instead of dilating the thresholded tag by 5 voxels,dilate it by 2 voxels, then subtract the thresholded tag and the dilatedpart from a 10×10×10 voxel windows. Since one only dilates 2 voxels, thetwo sides of the fold after the subtraction should still have a thinregion connected with each other. Region growing can find the connectedcomponent, an identifier is assigned to these connected voxels andmatting is performed again, with an added condition that there is nochange to the type of any voxel that has been identified as connectedtissue. In this way, one can avoid having true connected tissue beingreplaced by known air background. FIGS. 8( a)-(d) shows the result of atag on fold situation. The folds are the 2 linear structures 81, 82going from top-left to bottom-right in both the left side images, FIGS.8( a)-(b). The bright region 83 in the FIG. 8( a) on one of the folds istagged material. This material is removed after applying an algorithmaccording to an embodiment of the invention to obtain FIG. 8( b). Thehole 84 remaining on the tag can compliment the alpha foreground matthat was retrieved when polyp detection is performed. FIG. 8( c) shows azoomed in view that the fold 85 is still connected. The fold is notbroken because of air surrounding the fold matting inwards. FIG. 8( d)shows the material 86 that was removed from the FIG. 8( a) resulting inFIG. 8( b).

Small Tagged Materials with Intensity Lower than Threshold

Some tagged materials with size less than 3 mm×3 mm×3 mm have intensitylower than the foreground threshold. Their intensities are even lowerthan the true tissue intensity in some regions. Although these taggedmaterials usually do not effect detection of polyps bigger than 5 mm, amethod according to another embodiment of the invention can search forthese small tagged materials by creating sliding windows of 5 voxels×5voxels×5 voxels along the boundary between colon lumen and the colonicwall. If a voxel is a local maxima in that window, it can be treated asa small tagged material.

Speed

Since the task was confined to matting an unknown region, the region isjust a narrow-band around the tagged material with a thickness of about5 voxels. This allows for fast computation by looking at a much smallernumber of voxels than the whole 3D dataset.

Results Thin Cloth

An algorithm according to an embodiment of the invention can identify“fake” tissue which, due to partial voluming, has the same intensity asthe true soft-tissue, while retaining the real soft tissue near thetagged materials which might be misinterpreted as “fake” tissues. Thisis a classification task and the unknown region is treated as a locationwhere the voxels might be a mixing of either air with tagged materialsor tissue with tagged materials. FIGS. 9( a)-(d) shows the results ofclassification, foreground extraction, and background recovery. FIG. 9(a) is the original image, FIG. 9( b) shows the image with the foregroundextracted, FIG. 9( c) shows the image with the background recovered, andFIG. 9( d) shows the image with the background classified into tissueand air. The classification accurately classifies the background intotrue tissue and true air. The extracted foreground shows theprobabilistic map which indicates the blending fraction of foregroundand background. The background recovered looks natural. A simpleGaussian blurring can make the recovered boundary resemble the colonwall.

3D Matting

An algorithm according to an embodiment of the invention can incorporate3D information to analyze voxel type and voxel weighted intensity.Consider an area that looks like fake tissue due to partial volume canactually be real tissue being enhanced. FIGS. 10( a)-(d) illustrate how3D can capture more accurate information then 2D. FIGS. 10( b) and 10(d)suggest that the voxels 101, 102 at the arrow appears to be partialvolume between air and tag, but are in fact real tissue, as shown atanother viewpoint, FIGS. 10( a) and 10(c), which clearly show that thevoxels 101, 102 are within the soft tissue region. A 3D algorithmaccording to an embodiment of the invention can combine all the featuresavailable in 3D to aid the classification. In addition, some voxels nearthe tag are actually enhanced but fall below the foreground thresholdand cannot be visually identified by expert eyes, but the 3D spherecaptures the information by assigning fraction of blending in therespective slides.

FIG. 11 shows 2 adjacent slices. The top-left, FIG. 11( a), is one sliceand the bottom-left, FIG. 11( c) is the second slice. The top-right,FIG. 11( b), is the alpha value 111 detected in FIG. 11( a), and thebottom-right. FIG. 11( d), is the alpha value 112 detected in FIG. 11(c). In the first slice, the visible tag is below the threshold intensityso it is not detected. But with an algorithm according to an embodimentof the invention, the 3D sphere can pick the low intensity tag andassign alpha values. Looking at the next slice, the alpha value ispicked up correctly because there is a high intensity tag in this slice.Without the 3D sphere, the tag appearing in the first slice will bemissed.

More Results

FIGS. 12( a)-(c) shows a case where there is a hole between the tag andthe tissue. FIG. 12( a) is the original image, FIG. 12( b) shows region123 classified as background, and FIG. 12( c) shows the extractedforeground 124. The hole is the object 121 pointed to by the lowerarrow. The voxels in the hole are classified as air instead of tissue.It has a stand-alone tag which should be on the tip of a fold. Thesecond arrow 122 points to the stand alone tag. An algorithm accordingto an embodiment of the invention successfully removes the partialvolume around it and extracts the correct foreground. The classificationwould be incorrect would be if the hole was marked as tissue (whitecolor in FIG. 12( b)).

System Implementations

It is to be understood that embodiments of the present invention can beimplemented in various forms of hardware, software, firmware, specialpurpose processes, or a combination thereof. In one embodiment, thepresent invention can be implemented in software as an applicationprogram tangible embodied on a computer readable program storage device.The application program can be uploaded to, and executed by, a machinecomprising any suitable architecture.

FIG. 13 is a block diagram of an exemplary computer system forimplementing a method for detecting tagged material using alpha mattingin computed tomographic colonography, according to an embodiment of theinvention. Referring now to FIG. 13, a computer system 131 forimplementing the present invention can comprise, inter calia, a centralprocessing unit (CPU) 132, a memory 133 and an input/output (I/O)interface 134. The computer system 131 is generally coupled through theI/O interface 134 to a display 135 and various input devices 136 such asa mouse and a keyboard. The support circuits can include circuits suchas cache, power supplies, clock circuits, and a communication bus. Thememory 133 can include random access memory (RAM), read only memory(ROM), disk drive, tape drive, etc., or a combinations thereof. Thepresent invention can be implemented as a routine 137 that is stored inmemory 133 and executed by the CPU 132 to process the signal from thesignal source 138. As such, the computer system 131 is a general purposecomputer system that becomes a specific purpose computer system whenexecuting the routine 137 of the present invention.

The computer system 131 also includes an operating system and microinstruction code. The various processes and functions described hereincan either be part of the micro instruction code or part of theapplication program (or combination thereof) which is executed via theoperating system. In addition, various other peripheral devices can beconnected to the computer platform such as an additional data storagedevice and a printing device.

It is to be further understood that, because some of the constituentsystem components and method steps depicted in the accompanying figurescan be implemented in software, the actual connections between thesystems components (or the process steps) may differ depending upon themanner in which the present invention is programmed. Given the teachingsof the present invention provided herein, one of ordinary skill in therelated art will be able to contemplate these and similarimplementations or configurations of the present invention.

While the present invention has been described in detail with referenceto a preferred embodiment, those skilled in the art will appreciate thatvarious modifications and substitutions can be made thereto withoutdeparting from the spirit and scope of the invention as set forth in theappended claims.

1. A method for computer-aided object classification, soft segmentationand layer extraction in computed tomographic colonography, said methodcomprising the steps of: providing a contrast enhanced computedtomography (CT) digital image of the colon, said image comprising aplurality of intensities associated with a 3-dimensional (3D) grid ofvoxels; finding a foreground region of voxels with an intensity higherthan a pre-defined foreground threshold; creating a 3D trimap of saidcolon wherein said image is segmented into said foreground region, abackground region, and an unknown region between the foreground andbackground; starting from said background, extracting successive layersof said unknown region until said foreground region is reached, andclassifying each extracted layer as background or foreground; andgenerating a foreground matte, a background matte, and an alpha matte,wherein alpha indicates a mixing ration of foreground and backgroundvoxels.
 2. The method of claim 1, wherein said background region voxelsrepresent either regions of air in said colon, or regions of known softtissue in said colon.
 3. The method of claim 2, wherein creating a 3Dtrimap of said colon comprises morphologically dilating the 3Dforeground region one voxel at a time until background region voxels arereached, wherein a region through which said foreground was diluted is aregion of unknown type, and categorizing said background region voxelsas either air or soft tissue.
 4. The method of claim 1, whereinclassifying each extracted layer as background or foreground comprisesforming a 3D sphere around a candidate layer voxel and neighboring knownbackground voxels, and analyzing closest known background voxels todetermine whether said candidate layer voxel is a background air voxelor a background soft tissue voxel.
 5. The method of claim 4, furthercomprising calculating said mixing ratio alpha for each voxel in saidunknown region from alpha=(C−B)/(F−B), wherein C is the intensity ofsaid voxel, B is a background intensity, and F is a foregroundintensity.
 6. The method of claim 5, wherein analyzing closest knownbackground voxels comprises finding values of F and B that maximize asum of log-likelihoods L(C|B, F, α)+L(F)+L(B)+L(α) wherein${{L\left( {\left. C \middle| F \right.,B,\alpha} \right)} = {- \frac{{{C - {\alpha \; F} - {\left( {1 - \alpha} \right)B}}}^{2}}{2\sigma_{C}^{2}}}},{{L(F)} = \frac{{{F - \overset{\_}{F}}}^{2}}{2\sigma_{F}^{2}}},{and}$${{L(B)} = {- \frac{{{B - \overset{\_}{B}}}^{2}}{2\sigma_{B}^{2}}}},$wherein σ_(C) ^(r), σ_(F) ², σ_(B) ² are respective standard deviationsfor C, F, and B.
 7. The method of claim 1, wherein generating aforeground matte, a background matte, and an alpha matte, comprises:forming a 3D sphere about each candidate voxel in said region of unknowntype: obtaining from the background region within the sphere a number ofair voxels, a number of tissue voxels, their respective distance fromthe candidate voxel, and their variance; obtaining from the foregroundregion within the sphere a number of tagged voxels having an intensityabove a threshold for contrast enhancement, their respective distancefrom the candidate, and their variance; weighting the intensity of eachtissue voxel by its distance from said candidate voxel; weighting theintensity of each air voxel by its distance from said candidate voxel;weighting the intensity of each foreground voxel by its distance fromsaid candidate voxel; and calculating alpha fromC=αF+∫ _(Ω(air)) H _(air)(1−α(i))B _(air)(i)+∫_(Ω(tissue)) H_(tissue)(1−α(i))B _(tissue)(i), wherein C is the intensity of saidcandidate voxel, B_(air) is a background air intensity, B_(tissue) is abackground tissue intensity, F is a foreground intensity, Ω is a regionof integration in either air voxels or tissue voxels, i represents avoxel in the region of integration, and H is the Heaviside functionhaving a value of 1 inside the region of integration.
 8. The method ofclaim 7, wherein a foreground weighting function is:${F_{candidate} = {{\exp \left( {{- \beta_{1}}{\sum d_{i}}} \right)}{\exp \left( {{- \beta_{2}}/N} \right)}\frac{\sum{F_{i}/d_{i}}}{\sum{1/d_{i}}}}},$wherein F_(candidate) is the foreground weighted intensity at thecandidate voxel, F_(i) is the intensity value of the image at voxel iwithin the 3D sphere, d_(i) is the distance to the candidate voxel, β₁and β₂ are pre-defined parameters; a background weighting function is${B_{candidate} = {{\exp \left( {{- \beta_{1}}{\sum d_{i}}} \right)}{\exp \left( {{- \beta_{2}}/N} \right)}\frac{\sum{B_{i}/d_{i}}}{\sum{1/d_{i}}}}},$wherein B_(candidate) is the background weighted intensity at the voxellocation and B_(i) indicates the intensity value of the known backgroundwithin said 3D sphere.
 9. The method of claim 7, wherein a candidatevoxel is identified as a T-junction where the number of tissue voxels,air voxels, and a number of tagged voxels are greater than zero inside asphere of a predetermined small radius.
 10. The method of claim 7,wherein for tagged voxels laying on a haustral fold, the method furthercomprises: dilating said tagged voxels by a number of voxels less thanan original size of said morphological dilation; removing the taggedvoxels and dilated voxels from a voxel window whose side is greater thanthe original size of said morphological dilation; performing regiongrowing to find connected voxels; identifying the connected voxels; andperforming matting, wherein there is no type change of any voxel thathas been identified as a connected voxel.
 11. A method forcomputer-aided object classification, soft segmentation and layerextraction in computed tomographic colonography, said method comprisingthe steps of: providing a contrast enhanced computed tomography (CT)digital image of the colon, said image comprising a plurality ofintensities associated with a 3-dimensional (3D) grid of voxels;creating a 3D trimap of said colon wherein said image is segmented intoa foreground region, a background region, and an unknown region betweenthe foreground and background; forming a 3D sphere about a candidatevoxel in said region of unknown type; obtaining from the backgroundregion within the sphere a number of air voxels, a number of tissuevoxels, their respective distance from the candidate voxel, and theirvariance; obtaining from the foreground region within the sphere anumber of tagged voxels having an intensity above a threshold forcontrast enhancement, their respective distance from the candidate, andtheir variance; weighting the intensity of each tissue voxel by itsdistance from said candidate voxel; weighting the intensity of each airvoxel by its distance from said candidate voxel; weighting the intensityof each foreground voxel by its distance from said candidate voxel; andcalculating a mixing ratio of foreground and background voxels a fromC=αF+∫ _(Ω(air)) H _(air)(1−α(i))B _(air)(i)+∫_(Ω(tissue)) H_(tissue)(1−α(i))B _(tissue)(i), wherein C is the intensity of saidcandidate voxel, B_(air) is a background air intensity, B_(tissue) is abackground tissue intensity, F is a foreground intensity, Ω is a regionof integration in either air voxels or tissue voxels, i represents avoxel in the region of integration, and H is the Heaviside functionhaving a value of 1 inside the region of integration.
 12. The method ofclaim 11, wherein said foreground region comprises voxels with anintensity higher than a pre-defined foreground threshold for contrastenhanced voxels.
 13. The method of claim 11, further comprising,starting from said background, extracting successive layers of saidunknown region until said foreground region is reached, and classifyingeach extracted layer as background or foreground.
 14. A program storagedevice readable by a computer, tangibly embodying a program ofinstructions executable by the computer to perform the method steps forcomputer-aided object classification, soft segmentation and layerextraction in computed tomographic colonography, said method comprisingthe steps of: providing a contrast enhanced computed tomography ((CT)digital image of the colon, said image comprising a plurality ofintensities associated with a 3-dimensional (3D) grid of voxels; findinga foreground region of voxels with an intensity higher than apre-defined foreground threshold; creating a 3D trimap of said colonwherein said image is segmented into said foreground region, abackground region, and an unknown region between the foreground andbackground; starting from said background, extracting successive layersof said unknown region until said foreground region is reached, andclassifying each extracted layer as background or foreground; andgenerating a foreground matte, a background matte, and an alpha matte,wherein alpha indicates a mixing ration of foreground and backgroundvoxels.
 15. The computer readable program storage device of claim 13,wherein said background region voxels represent either regions of air insaid colon, or regions of known soft tissue in said colon.
 16. Thecomputer readable program storage device of claim 15, wherein creating a3D trimap of said colon comprises morphologically dilating the 3 Dforeground region one voxel at a time until background region voxels arereached, wherein a region through which said foreground was diluted is aregion of unknown type, and categorizing said background region voxelsas either air or soft tissue.
 17. The computer readable program storagedevice of claim 14, wherein classifying each extracted layer asbackground or foreground comprises forming a 3D sphere around acandidate layer voxel and neighboring known background voxels, andanalyzing closest known background voxels to determine whether saidcandidate layer voxel is a background air voxel or a background softtissue voxel.
 18. The computer readable program storage device of claim17, the method further comprising calculating said mixing ratio alphafor each voxel in said unknown region from alpha=(C−B)/(F−B) wherein Cis the intensity of said voxel, B is a background intensity, and F is aforeground intensity.
 19. The computer readable program storage deviceof claim 18, wherein analyzing closest known background voxels comprisesfinding values of F and B that maximize a sum of log-likelihoods L(C|B,F, α)+L(F)+L(B)+L(α), wherein${{L\left( {\left. C \middle| F \right.,B,\alpha} \right)} = {- \frac{{{C - {\alpha \; F} - {\left( {1 - \alpha} \right)B}}}^{2}}{2\sigma_{C}^{2}}}},{{L(F)} = \frac{{{F - \overset{\_}{F}}}^{2}}{2\sigma_{F}^{2}}},{and}$${{L(B)} = {- \frac{{{B - \overset{\_}{B}}}^{2}}{2\sigma_{B}^{2}}}},$wherein σ_(C) ^(F), σ_(F) ², σ_(B) ² are respective standard deviationsfor C, F, and B.
 20. The computer readable program storage device ofclaim 14, wherein generating a foreground matte, a background matte, andan alpha matte, comprises: forming a 3D sphere about each candidatevoxel in said region of unknown type; obtaining from the backgroundregion within the sphere a number of air voxels, a number of tissuevoxels, their respective distance from the candidate voxel, and theirvariance; obtaining from the foreground region within the sphere anumber of tagged voxels having an intensity above a threshold forcontrast enhancement, their respective distance from the candidate, andtheir variance; weighting the intensity of each tissue voxel by itsdistance from said candidate voxel; weighting the intensity of each airvoxel by its distance from said candidate voxel; weighting the intensityof each foreground voxel by its distance from said candidate voxel; andcalculating alpha fromC=αF+∫ _(Ω(air)) H _(air)(1−α(i)B _(air)(i)+∫_(Ω(tissue)) H_(tissue)(1−α(i))B _(tissue)(i), wherein C is the intensity of saidcandidate voxel, B_(air) is a background air intensity, B_(tissue) is abackground tissue intensity F is a foreground intensity, Ω is a regionof integration in either air voxels or tissue voxels, i represents avoxel in the region of integration, and H is the Heaviside functionhaving a value of 1 inside the region of integration.
 21. The computerreadable program storage device of claim 20, wherein a foregroundweighting function is:${F_{candidate} = {{\exp \left( {{- \beta_{1}}{\sum d_{i}}} \right)}{\exp \left( {{- \beta_{2}}/N} \right)}\frac{\sum{F_{i}/d_{i}}}{\sum{1/d_{i}}}}},$wherein F_(candidate) is the foreground weighted intensity at thecandidate voxel, F_(i) is the intensity value of the image at voxel iwithin the 3D sphere, d_(i) is the distance to the candidate voxel, β₁and β₂ are pre-defined parameters, a background weighting function is${B_{candidate} = {{\exp \left( {{- \beta_{1}}{\sum d_{i}}} \right)}{\exp \left( {{- \beta_{2}}/N} \right)}\frac{\sum{B_{i}/d_{i}}}{\sum{1/d_{i}}}}},$wherein B_(candidate) is the background weighted intensity at the voxellocation and B_(i) indicates the intensity value of the known backgroundwithin said 3D sphere.
 22. The computer readable program storage deviceof claim 20, wherein a candidate voxel is identified as a T-junctionwhere the number of tissue voxels, air voxels, and a number of taggedvoxels are greater than zero inside a sphere of a predetermined smallradius.
 23. The computer readable program storage device of claim 20,wherein for tagged voxels laying on a haustral fold, the method furthercomprises: dilating said tagged voxels by a number of voxels less thanan original size of said morphological dilation; removing the taggedvoxels and dilated voxels from a voxel window whose side is greater thanthe original size of said morphological dilation; performing regiongrowing to find connected voxels; identifying the connected voxels; andperforming matting, wherein there is no type change of any voxel thathas been identified as a connected voxel.